A study on correlation effects in two dimensional topological insulators 



(N 

o 

(N 
X> 

IX, 
IT) 



-f— > 

a 

o 
o 



m 
O 

(N 

m 

(N 
O 
(N 



Y. Tada, 1 R. Peters, 2 M. Oshikawa, 1 A. Koga, 3 N. Kawakami, 2 and S. Fujimoto 2 

institute for Solid State Physics, The University of Tokyo, Chiba 277-8581, Japan 
2 Department of Physics, Kyoto University, Kyoto 606-8502, Japan 
:1 Department of Physics, Tokyo Institute of Technology, Tokyo 152-8551, Japan 

We investigate correlation effects in two dimensional topological insulators (TI). In the first part, 
we discuss finite size effects for interacting systems of different sizes in a ribbon geometry. For large 
systems, there are two pairs of well separated massless modes on both edges. For these systems, 
we analyze the finite size effects using a standard bosonization approach. For small systems, where 
the edge states are massive Dirac fermions, we use the inhomogeneous dynamical mean field theory 
(DMFT) combined with iterative perturbation theory as an impurity solver to study interaction 
effects. We show that the finite size gap in the edge states is renormalized for weak interactions, 
which is consistent with a Fermi-liquid picture for small size TIs. In the second part, we investigate 
phase transitions in finite size TIs at zero temperature focusing on the effects of possible inter- 
edge Umklapp scattering for the edge states within the inhomogeneous DMFT using the numerical 
renormalization group. We show that correlation effects are effectively stronger near the edge sites 
because the coordination number is smaller than in the bulk. Therefore, the localization of the edge 
states around the edge sites, which is a fundamental property in TIs, is weakened for strong coupling 
strengths. However, we find no signs for "edge Mott insulating states" and the system stays in the 
topological insulating state, which is adiabatically connected to the non-interacting state, for all 
interaction strengths smaller than the critical value. Increasing the interaction further, a nearly 
homogeneous Mott insulating state is stabilized. 

PACS numbers: Valid PACS appear here 



I. INTRODUCTION 

Topological insulators (TI) and topological supercon- 
ductors (TSC) have been attracting great interest since 
the theoretical studies on the quantum spin Hall effect^— 
and their experimental realization^— They can be char- 
acterized by an energy gap in the bulk and gapless edge or 
surface states at their boundaries ! 10 ' 11 The direct relation 
between the bulk gap and the edge states, which is called 
the bulk-edge correspondence, is considered to be one of 
the most important properties in TIs and TSCs. These 
edge states are protected by topological properties of the 
bulk system, and are therefore robust against perturba- 
tions such as disorder or interactions which do not break 
the symmetries of the system . 12 ' 13 This suggests that, as 
long as the bulk remains gapful with non-trivial topology, 
the edge states exist at its boundaries. Correspondingly, 
in non-interacting systems, topological quantum phase 
transitions between topologically trivial and non-trivial 
phases require that the bulk gap is continuously closed. 
These gapped non-interacting systems which support the 
bulk-edge correspondence can be classified by ten distinct 
classes depending on their symmetries.— ~— 

However, it is not clear to what extent the bulk- 
edge correspondence holds in the presence of interac- 
tions. Interactions are expected to create novel corre- 
lated topological states such as TIs/TSCs without 
edge or surface states. For example, it has been proposed 
that strong correlations could realize a topologically non- 
trivial state without edge states when interactions be- 
come largeJ^ In this scenario, additional "zeros" in the 
Green's functions develop inside the bulk gap and finally 



merge with the poles which correspond to the edge states. 
Thus, the resulting state does not possess any edge states 
although it remains gapped in the bulk. Among the 
models exhibiting topological insulating states, the Kanc- 
Mele-Hubbard models have been extensively studied and 
the obtained phase diagrams show that the TI states ex- 
tend to a wide region.— ~— For large values of the Hub- 
bard interaction U, magnetic states related to the honey- 
comb lattice structure are stabilized. In between, a pos- 
sible realization of topological states without edge states 
has been proposed^ Correlation effects have also been 
investigated in one dimensional systems £ft2IL It has been 
shown that, in time-reversal symmetric TSCs with inter- 
chain interactions, the ground states are characterized 
by fermion parity in addition to the time-reversal sym- 
metry. These symmetries can distinguish different topo- 
logical states and the ground states arc labeled by topo- 
logical numbers in Zs instead of Z for non-interacting 
systems. Furthermore, it was pointed out that in some 
cases the edge states actually disappear as interactions 
are increased without any gap closing. 

Besides the possible realization of novel states in TIs 
and TSCs, strong interactions affect their quantitative 
properties. This is also the case for finite size systems, 
where the topological character would be smeared. Ac- 
tually, it was pointed out that the edge states can be 
gapped due to tunneling processes between the two sides 
of the sample, and physical properties strongly depend on 
external parameters such as temperature and magnetic 
fields.— ~— Finite size effects have been investigated also 
in the context of dimensional crossovers)^— However, 
finite size effects in the presence of interactions are still 
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not well understood, and properties of finite size TIs and 
TSCs would be sensitive to interactions in addition to 
external parameters. Moreover, correlation effects be- 
tween the two sides of the system could have significant 
scattering processes which might lead to novel inter-edge 
correlated states.— ~— 



In this article we discuss correlation effects in two di- 
mensional topological insulators in a ribbon geometry at 
zero temperature. In the first part, we discuss finite size 
effects in interacting TIs which generate a gap in the 
edge statesi22r— For large systems with L y ^> £ti, where 
L y is a width of the system and £ti is a characteristic 
localization length of the edge states, the massless edge 
states at each side are well defined and can be described 
by Tomonaga-Luttinger liquids (TLL)i22r— We study ef- 
fects of electron tunneling between the two edges by the 
standard bosonization method, and evaluate the finite 
size gap in spin and charge sectors. 




FIG. 1: (Color online) Ribbon geometry: square lattice with 
open boundaries in {/-direction. 



Next, we concentrate on relatively small systems, for 
which the edge states are massive Dirac fermions. For 
this purpose we use the inhomogeneous dynamical mean 
field theory (DMFT)^2r— with the iterated perturbation 
theory (IPT) as an impurity solver i 54 i 55 To be concrete, 
we analyze a Hamiltonian which has a topological non- 
trivial kinetic term as proposed for the HgTe/CdTe quan- 
tum wells and include an additional local Hubbard inter- 
action. We show that the finite size gap is renormalized 
by interactions, which can be naturally understood based 
on a Fermi-liquid picture. 

In the second part of this article, we discuss the 
paramagnetic Mott transition within the inhomoge- 
neous DMFT using the numerical renormalization group 
(NRG) 56 ' 57 as an impurity solver. We show that correla- 
tion effects are significantly stronger near the edge sites, 
because the coordination number is smaller at the edges 
than in the bulk, as has been discussed in correlated sys- 
tems without bulk gaps42 r 49 ' 58 Therefore, especially near 
the Mott transition, the localization of the edge states 
around the edge sites which is a fundamental property 
in TIs is weakened. However, possible inter-edge Umk- 
lapp scattering due to finite size effects does not seem to 
be important, and the "edge Mott insulating states" are 
not found. We show that even with strong correlations 
the ground state is a topological insulating state possess- 
ing edge states, which is adiabatically connected from 
the non-interacting system for all interaction strengths 
smaller than the critical value. When the interaction is 
stronger than the critical value, a nearly homogeneous 
Mott insulating state is stabilized. This is due to quan- 
tum tunneling of the electrons from the bulk to the edges, 
which was previously discussed in correlated metals. Fi- 
nally, in the appendix of the article we present details on 
the NRG calculations. 



II. FINITE SIZE EFFECTS 

In this section, we discuss finite size effects^r— in two 
dimensional TIs of different sizes in the ribbon geometry, 
see Fig. [T] When its width L y is large compared to 
the characteristic localization length of the edge states, 
L y ^> £rii the edge states, which are localized on the two 
sides, are well separated. Therefore, we can focus only on 
them and study correlation effects by using the standard 
bosonization technique. In this approach, forward scat- 
tering which has infrared singurarities is fully taken into 
account, while other interaction terms are perturbatively 
treated. We note that this consideration is based on the 
bulk-edge correspondence, which allows us to forget the 
details of the bulk of the system. 

For small system size where L y 3> £ti is not satisfied, 
the tunneling between the two sides is no longer a pertur- 
bation and the edge states are massive Dirac fermions. 
In this case, we take into account electron correlations 
by means of the inhomogeneous DMFT for a Hubbard 
model which will be introduced later. 

We first discuss relatively large size TIs in the next sec- 
tion, which is followed by the discussion of small systems 
in section Hi Bl 



A. Large systems: bosonization approach 

It is known that in the presence of interactions, each of 
the edge states, which are localized at the two sides, can 
be described in terms of the helical Tomonaga-Luttinger 
liquid (HTLL)i 12,13 However, when the two sides are cor- 
related, the edge states are usual Tomonaga-Luttinger 
liquids (TLL) with spin-charge separation.— ~— In this 
section, we investigate large size TIs by means of the 
standard bosonization approach ; 59 ' 60 focusing on inter- 
edge tunneling processes. 
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Fcrmion operators for the edge states, localized at the 
edges j = 1,2, are expressed as 

* jr<T (x) = e irkpX ij jr<T (x) 

for left(r = — or L) and right (r = + or R) movers with 
pscudo spins a =t)4-j and kp being the Fermi wavenum- 
ber of the edge states. The effective Hamiltonian for the 
two edge states is given by 

H = H 1 +H 2 + H ch + H f + H f , + H f - + H', 

Hi = -iv J dxlipl^dxTpiR-x - ip\ L ^d x i)iL\\-, 

H 2 = +iv J dx{ifjl L1 .d x ip2Lt ~ V'LyA'^flj.], 

H c h = 9ch J dx{p im p 1R f + PiliPili + (1 -» 2)], 

Hf = 9f J dx[p 1R ^pi Li + p2LtP2Rl], 
Hf = gf> J dx[pi R - t p 2 L^ + PlLlP2Rl\, 
H f"=9f" j dx[piRXP2R\ + PlLlP2Lt], 



where p, jra = ip^ r(J ij)j T(T and v is the Fermi velocity for 
the edge states. H lt 2 are kinetic terms for the edge states 
localized on edge-1 and edge-2 in the ribbon geometry, 
H c h ~ Hf are forward scattering terms, and H' includes 
other scattering processes. In the terminology of usual 
spin 1/2 fermions on a chain, H c h, Hf, Hf , and Hf" cor- 
respond to <?4|i, <?2_l 5 <72||, and g^i. terms, respectively. We 
note that, in the present study, SU(2) symmetry for spins 
is not required, and generally gf > gf>,gf>> holds. Al- 
though the indices R/L are redundant, we keep them 
so that the relation between the usual spin 1/2 fermions 
and the present system is clear. We also note that the 
parameters in the above Hamiltonian should be regarded 
as those including correlation effects in the bulk system 
from which the effective one dimensional model is de- 
duced. 



Possible scattering processes in H' are 



H = H\ + My + Hhf + H s f + H e f + H u + H u > , 

H x = \J dx[<pl R ^2L^- l2k ^ + 4 Ri ^Ll e ~ l2kFX + ( h - c -)L 

H\> = A' J dx[lp\ m Tp2Ri ~ ^iL^Lf + (h.c.)], 

H h f=g hf J ^[^^^^^^^^^-^.^^^t^^^^^^+^fe^^ 

Hef =9ef J dxllpl^lpl^tpiL^Rf + V^RfV^^i^LtL 
Hsf=gsf J dx^l^l^iL^Ri + ^Ih^Il^R^L^ 

H u = g u J dx[i>\ m {x)ilj\ m {x + a )ip 1Li {x + ao)*Pi H (x)e-' l4kFX + ^\ h ^\ L ^2R^2R^ kFX + (h.c.)], 
H u , = g u , J dx[^\ R1; {x)^\ m {x + a )i(j2Lt(x + a )^ 2 Lt(x)e^ likFX + tpl Rl ^l Ri ^iLi^iue- l4kFX + (h.c.)], 

I 



where ao is the lattice constant. H\ is a spin-conserving 
tunneling term, and H\> is a spin-flip tunneling term 
which can be nonzero if the bulk Hamiltonian does not 
have the spin rotation symmetry. Hhf,H e f, and H s f 
describe " hclicity-flip" , "edge-flip" and "spin- flip" pro- 
cesses, respectively. H u and H u > arc Umklapp terms. 
Hhf,H e f and H u i correspond to <?3j_,gij_ and g 3 j[ terms 
in the spin 1/2 fermions. Obviously, they are invari- 
ant under the time- reversal transformation tpj R /L^(x) — > 
TpjL/Riix) and ipj R/Li {x) -> -ip jL/m (x). 



To investigate the effects of H', we use the standard 
bosonization approach i 59 ' 60 Bosonization relations are 



U 



7 e -ir<t>ir 



4 



where Uj ra are the Klein factors. We introduce 

<t>c = (4>lLi + 4>lRt + 4>2Lt + 4>2R\)/2\2, 
Oc = {4>lLl - 4>lRt + 4>2L-t ~ 02_R;)/2v / 2, 

4>s = (-(f>iLi + 4>iR\ + ^2Lt _ 4>2Ri)/2v2, 

6s = (-<filLl - <f>lBt + 02Lt + ^2Ki)/2v / 2. 

These variables satisfy 

[<f>» (x), d6 v {y)] = i-K5 llv 6{x - y), 

and the canonical momentum is 11^ = (l/^dd^. 

Then, the kinetic term in the Hamiltonian Hk = H\ + 
H 2 + H c h + Hf + Hfi + Hfn can be diagonalized and the 
bosonized Hamiltonian is, 



with g c h = gch/K and g = g/2-K for the other terms. H c 
and H s describe collective charge and spin excitations, 
respectively. As mentioned above, because gf > gti, K a 
is not fixed to unity and generally K s > 1. If all the 
inter-edge correlated terms are neglected, H c and H s are 
essentially identical, and the edge states are nothing but 
two copies of the HTLL localized at each edge. We note 
that, if the bulk system is simply renormalized by a renor- 
malization factor z in the presence of interactions, v is 
proportional to z and g is to z 2 . Therefore, deviation of 
K SyC from unity would be suppressed by a factor z. 



where 



H c 



H c + H s , 





f dx 


- 1 


2tt j 


-K~c 






- 1 




f dx 




2tt j 







Kc(d9 c ) 2 
K s (d9 s f 



[(v + g C h + gf - 
x (v + g ch + g r 

v + g c h+ 9f" - 



■9f+9f) 
9f ~ 9f> 



v + 9ch + 9f" - 
[(v + g c h - 9f - 
x (v + g ch - gf» 

v + 9ch - 9f" - 



1/2 



9f + 9f 
■9f-9f>) 
-9f+9f')] 1/2 , 

- 1 1/2 

' 9f ~ 9f> 1 



9ch - gf -gf + gf 



H' is also bosonized as 



Hx- 

H h f = - 
H ef = u 
H s f = - 
H u 
H,.i 



iX f iX 

Ui R ^U 2 l\ I dx sin(x/2((f> c + <f> s ) - 2k F x) H U 2 r\Uxli / dx sin(v / 2(</) c - <j> 8 ) - 2k F x), 

7ro J na 



— U 2L \U 1Li f dxsmV2(<j> a +6 s ) + —U im U 2 Ri I dxsin V2(^ s - 6 S ), 

— u i — / dx cos( V&4> c — ik F x), 
{2ira y J 



{ dxcos\/8<f> s 
liraoy J 



(2 

2g^/ 

M (27ra ) 2 
2.9„ 



(27ra ) 5 



dx cos 



dx[cos(y/8((/> c + 9 S ) - 4k F x) + cos(v / 8(</> c - 6 S ) - ik F x)], 



^a )2 J dx l cos (^(<l>c + 4>.) ~ 4k F x) + cos(V8(0 c - (f>„) - 4Js F x)) f 



where u = U2rxU\r^U\l±U2l-x ■ Hx and Umklapp terms mix the charge and spin sectors. 



5 



The tunneling Hamiltonian H\ includes 2kp- 
oscillating factors, because it describes hybridizations 
between edge states on the two sides with the same 
wavenumber k. In the following, we focus on a special 
filling, kp = 0, to evaluate the finite-size gap induced by 
the tunneling Hamiltonian. When kp = 0, the chemical 
potential is exactly at the Dirac point of Hi ^ and 
all the oscillating factors in the Hamiltonian become 
unity e l2kFX = e l4kpX = 1. Then, we can evaluate the 
perturbations by looking at their scaling dimensions 
for H,j., 



D x 


= (K c + 


K s )/2, 


Dy 


= {K S + 


l/K s )/2, 




= 2K C1 




Def 


= 2K S , 




Dsf 


= 2/K s , 




D u 


= 2{K C 


+■ 1/K.), 


D u . 


= 2{K C 


f K s ). 



From these equations, one sees that D\ is the small- 
est scaling dimension around the non-interacting point 
K c = K s = 1, and therefore, the spin-conserving tunnel- 
ing H\ would be the most relevant perturbation in H' . 
Therefore, the edge states become gapped by H', and 
they are no longer the TLL in finite size systems when 
kp ~ 0. We note that, since the tunneling terms H\ 
and H\i have non-zero conformal spins, they contribute 
to renormalization group flows of f , H e f and H s f , re- 
sulting in complicated flows^r— However, because H\ is 
the most relevant perturbation, we focus only on H\ and, 
in the following, simply neglect all the other perturbation 
terms in H' . 

Under this approximation, we can estimate the finite 
size gap by using the standard analysis^ When A is neg- 
ative and large enough, conditions \[2 (<f> c ± <f> s ) — tt/2 
give the minimal energy and the Hamiltonian would be- 
come 



A 

2 

ira 

2A 
(ira ) 



A 
na 



2{4> c -4>sY 



<4> 2 c + 4> 2 ), 



where A = 27raoA. Here, we have used a representation of 
the Klein factors where C/iflfL^Lt = U2r\U\l\, = —i can 
be satisfied, like U\r\ = u\ g) oi, U\li = 1 <S> 02, U2R^ = 
0"3<8>cri, Uil-\ = <J2®<J\, where {<Ji\ are the Pauli matrices. 
4> c has been shifted by n/2y/2. Then, the spin and charge 
sectors are decoupled and the gaps induced by H\ arc 
respectively evaluated as 



i/(2-£>a) 



w 



where w is an effective band width of the edge states 
which would correspond to the bulk gap Ati- A SjC are 



enhanced by the forward scattering. We stress that, even 
for the weak interactions K s c ~ 1, both of the spin and 
charge sectors are gapped by the tunneling H\. There- 
fore, the TLL description fails for the low energy edge 
states in the finite size TIs when kp ~ 0. Neverthe- 
less, correlation effects peculiar to one-dimensional sys- 
tems appear as a difference between the charge gap and 
the spin gap, ^ A^, which characterizes the spin- 
charge separation. We stress that the charge gap can be 
much smaller than the spin gap in strongly correlated 
regions provided that the edge states can be treated as 
one-dimensional systems. It should be an interesting is- 
sue to experimentally detect the spin-charge separation 
raised by the one-dimensional correlation effect from the 
measurement of the Hall effect and the spin Hall effect 
which are governed by the edge states. 

In the above derivation of A s c , we have simply ne- 
glected other interaction terms in H' . However, they can 
be important in a strongly correlated regime where K s c 
is away from unity. As mentioned above, Hy , H s f , H u 
describe spin flip processes which are difficult to take 
into accout in numerical calculations. On the other hand, 
the inter-egdge interactions Hhf, H u > are spin-conserving 
Umklapp terms and might lead to the possible "edge 
Mott insulating states". The effects of the inter-edge 
Unklapp scattering are discussed in Sec. Mil 

At the end, we note that all the parameters in the ex- 
pressions for the gaps, K s ^ c , A and w, are derived from the 
two dimensional bulk system and would be renormalizcd 
by the interactions in the bulk. This renormalization 
cannot be calculated from the one dimensional model, 
and rather should be included in the initial values of the 
parameters. 

In the next section, we analyze small TIs by an alter- 
native approach which is based on a Fermi-liquid picture 
for the TIs. 



B. Small systems: DMFT+IPT calculations 

When the system size is small enough so that the tun- 
neling between the edges is no longer a small perturba- 
tion, analyses based on well-defined masslcss modes are 
not applicable anymore. In this section, we consider cor- 
relation effects in small size systems whose edge states 
are massive Dirac fermions. 

To be concrete, we consider the following Hamilto- 
nian defined on a square lattice in the ribbon geometry 
(Fig- ID, 



H = H 



BHZ 



H 



hit 



(1) 



where £/bhz is topologically equivalent to the Hamil- 
tonian for the HgTe/CdTe quantum wells^^^i H int 
is an interaction term between the electrons. The 
Hamiltonian contains two orbitals described by Cj = 
(Ciif, Ci2f, c iU ,Ci2i) T where the indices 1, 2(t, j.) label or- 
bitals (pseudo spins). The kinetic energy can be written 
as 
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H B HZ - ^2 CjHijCj, 



Hi. 



Hi 



Hij 
H*j 



t'[i(S. 



MoSij — t(5ij±x + Si,j±y) 



i,j+x 



) + Si, j - 



t'[i{5, 



5i,j-x) + $i,j+y 

k,j+y] -M Q 5ij + t(SiJ± X + SiJ±y) 



(2a) 

(2b) 
(2c) 



The periodic boundary condition is imposed for the a;- 
direction, while the system has a finite width L y for the 
2/-direction. As will be discussed later, this Hamiltonian 
is particle-hole symmetric when the filling is ri\ + n-2 = 1 
per spin-direction. Regarding the interaction term, we 
consider the simplest Hamiltonian 



H, 



(3) 



where nn a = c ila cu a . The total Hamiltonian is a Hub- 
bard model with inter-site inter-orbital hopping t 1 corre- 
sponding to the spin-orbit (SO) interaction in the present 
article. In the following, we take t as the energy unit 
t = 1, and the SO interaction t' is fixed to be t' = 0.25. 
Furthermore, we set the filling to Tin + mi = 1 per 
spin-direction for all sites. We also fix L y = 10 and 
Mo = —1.0, if not specified, for which the finite size gap 
is Af ~ O.IAti in the absence of the interaction. As 
discussed in the previous section, in finite size TIs, in- 
frared singularities are cut off and descriptions based on 
the TLL fail to capture gapped edge states when their 
Fermi wavenumber is small. Instead, a Fermi-liquid pic- 
ture of correlated band insulators would give a correct 
understanding where the edge states are massive Dirac 
fermions, as will be discussed in the following sections. 



1. Non-interacting case 

Before discussing correlation effects, we want to show 
some basic properties of the non-interacting system. In 
Fig. [21 the dispersion relation for M = —1.0 calculated 
from i?BHZ is shown. The dispersion is particle-hole sym- 
metric and edge states exist for < \Mq\ < 4. The bulk 
gap is Ati — 1-0 corresponding to t' = 0.25. The edge 
states have a small but finite gap as discussed in previ- 
ous studiespSr— The finite-size gap Af depends on the 
parameter Mq and it is oscillating when Mq is tuned as 
shown in Fig. [3] for L y = 10, 20 and 40. The number 
of the oscillation for < \Mq\ < 4 is 2(L y — 1), which 
is related to the number of the eigenstates of Hbhz for 
a fixed k x . This oscillation is specific to lattice systems, 
and we have to be careful about this when discussing cor- 
relation effects in the following sections where oscillations 
in Af can also be seen. Oscillations in bulk gaps with re- 
spect to system sizes have been discussed in many articles 




FIG. 2: (Color online) The dispersion relation for Mo = —1.0. 
The edge states spectrum is enlarged in the lower panel so that 
the finite size gap can be seen. 



and oscillatory dimensional crossovers are seen in some 
classes of TIs and TSCs.— _ — However, in our study, the 
non-interacting system remains topologically non-trivial 
as long as < \Mo\ < 4. Af is exponentially decreasing 
as L y is increased, and its size dependence is roughly 



A f ~ Af exp (-L y /£, 



4°> ~ t'/M , 
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FIG. 3: (Color online) The finite size gap as a function of Mo 
for L y = 10, 20, 40 from the top to the bottom. For < Mo < 
4, Af depends on Mo in the same way as for —4 < Afo < 0. 



where characterizes the localization of the edge states 
around the edge sites in the non-interacting system. Afo 
is defined by a combination of the parameters in -ffeHZ- 
Exact expressions for these quantities in the continuum 
limit can be found in Zhou et al.— . 

Naive expectation based on the Fermi-liquid picture 
together with an assumption that in the presence of the 
interaction U the low energy states are renormalized by 
a single renormalization factor z, results in a finite-size 
gap 



A f0 (C7) ~ zA l0 (U = 0), 
ZTi{U)~zi!/zMo~i 



(0) 
TI ' 



(4a) 
(4b) 



Renormalizations of if and Mo would cancel out and 
£ti{U) would be equal to . Therefore, our simplest 
expectation is that the amplitude of the finite size gap 
is renormalized by z while £ti is not, although the legit- 
imacy of the above assumption is not so clear. In the 
following sections, we investigate correlation effects and 
show that the naive expectation actually holds in the 
presence of the Hubbard interaction within the DMFT. 



2. Inhomogeneous DMFT+IPT 

In this section, we discuss the correlation effects by 
means of the inhomogeneous DMFT+IPT. The inho- 
mogeneous DMFT is an extension of the DMFT— to 
inhomogeneous systems and has been applied to, e.g. 
metallic surfaces ; 44 ' 45 heterostructures^^— and optical 
lattices.— ~— We note that, in the present study, the edge 
states for small L y arc massive and not the massless 
HTLLs as discussed in the previous sections. The whole 
two dimensional system including the edge degrees of 
freedom is supposed to be a gapped Fermi-liquid. There- 
fore, we expect that inhomogeneous DMFT can capture 
the main properties of the edge states in the presence of 



interactions, although it is known that DMFT fails to de- 
scribe the low energy TLL behavior in one dimension.— 
To calculate the selfcnergy, we solve the following self 
consistent equations for y = 1, 2, ■ ■ • , L y , 



— 'Y^G(u,k x ,y,y) 



+ E 00 /(a;,2/), 



where a = (la), Q is the cavity Green's function and L x is 
the number of k x points, which is taken to be sufficiently 
large. The lattice Green's function G aa i(uj,k x ,y,y') is a 
AL y x AL y matrix, where the system size is fixed as L y = 
10 in the present study. We note that G is spin-diagonal 

and J2k x G W°( u ' k x>y>y) = Sfc„ Gi<7,i'a(u,k x ,y,y) is 
satisfied in our ribbon geometry, and therefore we often 
suppress the spin indices throughout this article. 

The selfcnergy of the effective impurity problem is cal- 
culated by means of the IPT which can compute the self- 
energy as a function of real frequency a; . 54 ' 55 Within the 
IPT, the selfenergy is evaluated in the following way: The 
first order contribution is S^; /(r , (y) = Uni a (y)Sw5 crcr ' 
where ni a (y) is the electron density for orbital I and spin 
a at site y. We note that n i<y(v) = 1 f° r au the sites. 
The second order is 
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where a, b,c,d = (la) and T is an antisymmetric bare ver- 
tex. / is the Fermi distribution function. In the present 
model, r^'j^ = U,Tfofo = —U, and all the other el- 
ements are zero. The Green's function Q is given by 
= Giik'o' ~ UniaSu>Saa>- The selfenergy is inter- 
polated from the weak coupling to the atomic and high 
energy limit in the IPT calculations, 



,(i) 



(2) 



i 



where A w = 5 w ni(l - ?ij/2)/nj(l - nj/2) + (1 - 6 W ) 7 
Bui = Su'2(l — ni)/ni(l — hi/2), and hi is calculated 
from Q . We use a simplified scheme which has been suc- 
cessfully applied to correlated metals away from the half 
filling, in which the effective chemical potential is simply 
determined from E(w = 0) = E'^O)^ 

In Figs. [2]arid[SJ we show the selfenergy and the local 
Green's functions for half of the system width for [7 = 4. 

As will be discussed later, the local selfenergies for the 
other half are ImEn .22(^7 L y — y + 1) = ImEix 22(^5 y) 
and ImEi2,2i(w,£y - V + 1) = — ImEi2,2i(w, y). The 
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FIG. 4: (Color online) Selfenergy for U = 4 and L y = 10. 
Red, green, blue curves are En, £12,21 and E22, respectively. 



local Green's functions show the same symmetry prop- 
erties. We see that, the diagonal elements of the selfen- 
ergy are similar to those of simple band insulators; e.g. 
Im£ ii^Cw) seem to have gap like structures. They be- 
have like ImSn j 22(w, y) ~ — b(y)ui 2 for small w, where the 
prefactors b(y) are much smaller than those in metallic 
systems. This is because the density of states has only 
small weight around uj = which corresponds to the edge 
states. On the other hand, although ImGi2,2i is compa- 
rable to ImGn^ in magnitude at least for y = 1 around 
(jj r*j Ati, £12,21 is much smaller than £11,22- In Fig. [6l 
we also show the energy spectrum for U = 4 defined by 
A(k x ,uj) = -(l/L y )J2ylm[G 11 (k x ,y,u) + G 2 2(k x ,y,uj)}. 
Compared to Fig. ® the edge states around uj = are 
renormalized and high energy structures are smeared. 
We note that, although the massive Dirac states are topo- 
logically protected, namely they are free from the single 
particle backscattcring and the intra-edge Umklapp scat- 
tering in the spin conserved systems, they have a long 
but finite lifetime due to the forward scattering. A di- 
rect comparison with a non-perturbative impurity solver, 
the NRG, is shown in Appendix I A 31 

The finite size gap Af as a function of U can be calcu- 
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FIG. 5: (Color online) The local Green's functions for y = 
1 ~ 5 from the top to the bottom when (7 = 4 and L y — 10. 
The red, green, and blue curves are for G11, Gi2,2i, and G22, 
respectively. 
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FIG. 6: (Color online) The energy spectrum A(k x ,u)) for U = 
4 and L v = 10. 



latcd from an effective Hamiltonian 

H(uj) = H + Z(lj). 

For |w| < Ati, Im£(o;) is much smaller than that for 
t' = and can be neglected when evaluating the finite 
size gap. We can identify the positions of zeros in u> — 
H(ui), which correspond to the energy eigenvalues of the 
edge states. As shown in the upper panel of Fig. [71 Af 
is an oscillating function of U. This oscillation has its 
origin in the non-interacting system in Fig. [31 because 
constant terms in Rc£n.22 play the same role as Mq in 
the inverse of the Green's functions. The bulk gap Ati is 
simply estimated from the next lowest excitation energy 
of H at k x = 0, and Z = (1/ L y )J2 y z (y) is also shown 
in the figure. We note that, compared to Ati (£7), the 
renormalization of Af is rather moderate, because z(y), 
which is defined by the derivative of £ at uj = 0, is larger 
than that at w ~ Ati- To extract the correlation effects 
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FIG. 7: (Color online) (Upper panel) The finite size gap Af , 
the bulk gap Ati, and the y-averaged renormalization factor 
Z as functions of U for L y = 10. (Lower panel) The ratio of 
the finite size gap 5f and Z as functions of M. 
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FIG. 8: The elementary excitations in two dimensional TIs 
for different system sizes when the Coulomb interactions are 
weak. 



longer applicable. 

In this section, we have discussed the correlation ef- 
fects for small size systems, while the previous section 
was dedicated to rather large systems. For these two lim- 
iting situations, the physical pictures arc very different; 
the edge states are renormalized massive Dirac fermions 
for small systems, instead of the gapped spin and charge 
collective modes for large systems. However, it would 
be possible to smoothly connect them to each other as 
the system size L y is tuned. If one starts from a small 
size where the single particle spectrum has a gap, one 
should calculate the susceptibilities to see the gaps in 
the spin and charge collective excitations. These gaps in 
the susceptibilities should be evaluated with the vertex 
corrections being taken into account. Although such cal- 
culations are beyond the present DMFT calculations, a 
crossover would take place at some length scale of the sys- 
tem. Such crossovers between two or three dimensional 
physics and one dimensional physics with respect to en- 
ergy scales are also seen in ladder systems.— In Fig. 
we summarize the elementary excitations in two dimen- 
sional TIs for different system sizes when the Coulomb 
interactions are weak. In the limit L y — > oo where the 
inter-edge correlation is negligible, the "gapless collective 
excitations" correspond to the HTLL. 



from this oscillating Af, we replot it as a function of a 
parameter M corresponding to M . We define 

M(U) = Mo + ^-^ReEn^O,?/) 

y y 

-ReS 22 (w = 0,y)], 
St(M) = A f (M(U))/A f (M ;U = 0)\ Mo=M , 

where VWo = —1.0. In the definition of i5f , we have identi- 
fied the non- interacting parameter Mq as M. . In the lower 
panel in Fig. [Jj we show S{(A4) and Z(A4). We see that 
5[(A4) and Z(Ai) well coincide, which means that our 
naive expectation Eqs. Paj) and (|4b[) based on the Fermi- 
liquid picture actually holds within the DMFT+IPT cal- 
culations. We note that, for weak interactions, site de- 
pendence in is not strong and Af can be well un- 
derstood in terms of their averages, M. and Z. If the site 
dependence in £ is significant, this interpretation is no 



III. MOTT TRANSITION 

In this section, we discuss the paramagnetic Mott 
transition for a topologically nontrivial system by us- 
ing the inhomogencous DMFT. We use the same model 
Hamiltonian as in the previous section, Eq. (fl]), and 
suppress possible magnetic instabilities. In previous 
studiesi£ii&2&2&~— it was discussed that for TIs and 
TSCs with strong interactions novel ground states might 
be realized, in which the bulk remains gapped with a 
non-trivial topological number while the edge states are 
suppressed by interaction effects. Such states are called 
"edge Mott insulating states", and can be explicitly de- 
fined in terms of the Green's functions^ Let us briefly 
consider a semi-infinite TI, which is built up as a het- 
erostructure of the vacuum and the TI having one bound- 
ary. In the absence of interactions, poles are supposed to 
be in the Green's functions corresponding to the edge 
states, while in the presence of interactions there might 
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be additional zeros as was discussed by Gurarie— . In 
the latter case, the edge states are suppressed and the 
the ground state can be characterized as an "edge Mott 
insulating state" . Although this discussion cannot be di- 
rectly applied to finite size systems with two boundaries, 
the selfenergy must have poles inside the bulk gap to re- 
alize the edge Mott insulating state, which characterizes 
the TIs with Mott insulating edge states. In this sce- 
nario, the finite size gap would be replaced by the Mott 
gap as interactions are increased. We note that, to obtain 
the edge Mott insulating states, the Umklapp scattering 
for the edge states is essential. Although, the Umklapp 
scattering of the edges at the same side is absent in our 
model, there are possible Umklapp scattering processes 
which involve both sides. In the edge state Hamiltonians 
introduced in Sec. Ill Al these processes are given by Hhf 
and H u >, while H u is absent. Such scattering processes 
might drive the edge states into a Mott insulating state, 
keeping the bulk a topological insulator. We note that, 
while H u i is characterized by scattering at different sites, 
Hhf describes scattering at a single site. Because there 
is some overlap of the edge state wave functions which 
are localized mainly at opposite sides, effects of Hhf arc 
taken into account within the inhomogencous DMFT. In 
addition, when Af/Axi is not so small and U is large 
compared to Ati, the inter-edge Umklapp scattering is 
expected to be important. 



A. Symmetries in Green's function and selfenergy 

Before discussing results obtained with inhomogeneous 
DMFT, we discuss some identity relations in G and £ 
arising from the particle-hole symmetry in the Hamilto- 
nian Eq. (fTJ). These identities are quite useful for check- 
ing numerical results and simplifying calculations, thus 
improving the accuracy. 

The particle-hole transformation, defined by a unitary 
operator C, is given as 

C-ila c i2cr' 
J 



Ci2a 



The condition that the Hamiltonian should be invari- 
ant under this transformation determines the chem- 
ical potential. The interaction in the model is 
given by H int = U Y^a nu^nu^. Thus, the chemi- 
cal potential should be /i = U/2 so that particle- 
hole symmetry is conserved. We note that our 
model belongs to the class DIII-Hamiltonians which 
describe the p ± ip superconductors such as the 

Rashba superconductors^ if we identify our particle- 
hole symmetry C as emergent particle-hole symmetry 
in superconductors^ This correspondence would allow 
us to interpret the off-diagonal selfenergy £12,21 as an 
anomalous selfenergy of the superconductivity if the 
emergent particle-hole symmetry is kept in the presence 
of interactions. 



From the particle hole symmetry we can directly derive 
identities in G(w) and E(w). The local retarded Green's 
function for site i satisfies 



Gr 2 ( 

@22 (' 



U,l) 



-G&(- 
-G&(- 



■U),l) 



(5) 



We note that the inverse of the matrix G(ui) = 
[Gui{uj,i, j)] also satisfies the corresponding relations, 
since G(w)- 1 • G(u) = G{u) ■ G(w)- 1 = 1. 

Thus, the selfenergy for a fixed filling fulfills similar 
relations, which arc directly derived from 



S fl (/i;w,*) 



[G^;w)\ 



= [[G*( W oj)}^ - [G"( Mo ; w)]^J + (/x - Moj 

where each of the non-interacting and interacting lo- 
cal Green's functions, the site-diagonal elements of the 
matrices G r (hq;uj) and G R (fi;uj), satisfies Eq.©. The 
chemical potential is /iq = and [i = U/2. Therefore, 
the local selfenergy should satisfy 

£fi(-w,£)-J" ^(-uj) 
£fi(-w,«) Ef 2 (-a;,i)-)ti 



(6) 



We see that the real part of £12,21 (w) is an odd function 
of ui and the imaginary part is an even function, and es- 
pecially Re£i2,2i(w = 0) = 0. Looking at the inverse of 
the Green's function G -1 = [ui — H — £] , the off-diagonal 
selfenergy £12,21 can be interpreted as a local hy- 
bridization between the orbitals, which does not add to 
the inter-site inter-orbital hopping t' in the Hamiltonian 
H . In this interpretation, £12,21 (w, y) are regarded as the 
induced " s-wave" hybridization between the orbitals aris- 
ing from the "p-wave" hybridization t' in the presence of 
the on-site interaction^ The odd frequency dependence 
in Re£i2, 21(^0 is related to the odd parity character of 
the inter-orbital hopping t' . 

The odd parity £12,21 in our system is in contrast 
to those in trivial band insulators, where inter-orbital 
hopping is local, s-wavc-like. It can be shown for this 
model that off-diagonal elements behave like £12 (a;, i) = 
Eia(— Therefore, Re£i 2 ,2i(w = 0, i) ^ can 
add to the local hybridization, which suggests that the 
bulk gap grows continuously as U is increased. This has 
been shown in the previous studies for topologically triv- 
ial correlated band insulators In the context of su- 
perconductivity, this would correspond to the BCS-BEC 
crossover for the s-wave pairing states i 70 ' 71 

We also note that this argument can easily be gener- 
alized to the non-local selfenergy, and the off-diagonal 
elements in the present model satisfy £12(0;, = 
— £12 (— oj, j, i)* for sites At u = 0, £12(0, i,j) has 
the same symmetry as the inter-site inter-orbital hop- 
ping and enhances it. We again stress that the odd- 
parity symmetries in the off-diagonal selfenergy are char- 
acteristic of our model. However, constant terms in 
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ReSi2(a;, i, j) do not appear in the lowest order in U, 
and ReS 12 (0, i,j) comes only from higher order terms in 
U. Although ReEi2(0, i,j) is expected to be small in the 
present model, sufficiently large ReEi2(0, might lead 
to a crossover with continuous growth of the bulk gap in 
the presence of on-site interactions. However, previous 
studies on the Kanc-Mclc Hubbard model taking into ac- 
count spatial correlations show a rcnormalization of the 
bulk gap^i, which implies that in our model Ei 2 (0, 
for i ^ j is not dominant, too. 

Other relations in G and E corresponding to an in- 
version transformation along the y-axis also exist. The 
inversion transformation I is defined by 



orbital-1 



where i = (—x,L y + 1 — y). Under this transforma- 
tion, the Hamiltonian satisfies XH{t - . ' , t i ■ ' )X _1 = 

11 22 12 21 

Hfaj' ,—tij' )■ Because Gn. 22 are even with respect 
to £12,21 while Gi2,2i are odd, the following relations hold: 



G?Mi) G? a (w,i) 
Sf!(w,i) Sg(w,i) 



Gf^w.i) -Gg(«,i) 

-G£(w,i) G&(w,i) 

EfiCw.i) -Ef a (w,i) 

-Ef 1 (w,i) Ef 2 (w,£) 



(7) 



(8) 



We note that, although Gi 2i2 i(u;,y) are comparable to 
G\i,22{u,y) in magnitude around the edges, they are 
strongly suppressed in the bulk. They become exactly 
zero at the center of the system if L y is odd, as can be 
seen in the above identities. 



B. DMFT + NRG calculations 

We will now discuss the correlation effects in small size 
TIs from moderate to strong coupling within the inho- 
mogeneous DMFT using the numerical renormalization 
group (NRG) 56 ' 57 as an impurity solver. The system size 
is fixed as L y = 10. If nothing else is stated, we use a dis- 
cretization parameter A = 3.0 and Nk = 5000 states are 
kept in the calculations. The effective impurity Anderson 
model for a single lattice site reads 

H\am{v) = H c (y) + H hyb (y) + H imp (y), 

Hc = X ala( ms )9u'<j( ms ) a i'<j{nis), 

ms,ll' a 

H hy b = X \^f} IJ hu'a{'ms)av cr (m,s) 

msAl' a 



Ha, 



+a\, a (ms)h* ul(T (ms)fi 

X! E LfLfl° + H tn t[f} a , fla] 




where the parameters are determined self-consistently 
within the DMFT loop and depend on the y-index of 



impurity 



orbital-2 

FIG. 9: (Color online) A schematic picture of the NRG chain 
Hamiltonian of the effective impurity Anderson model. 



the lattice. Details for the derivation of the parameters 
are given in Appendix IA II and IA 21 The Hamiltonian 
can be rewritten as a linear chain as shown in Fig. [5] To 
examine the effects induced by strong correlations, the 
NRG is expected to be a very accurate impurity solver, 
especially near the Fermi energy 56 ' 57 However, due to 
the SO interaction and the resulting inter-orbital hop- 
ping in the chain Hamiltonian, truncation effects can be- 
come large in these calculations. Therefore, the NRG 
must also be considered as an approximation even around 
the Fermi energy. (A discussion on truncation effects is 
given in Appendix IA 31 ) Furthermore, instead of using 
the usual logarithmic broadening of the delta-peaks in 
the Green's functions calculated by the NRG^ 7 - we use a 
Lorentzian broadening. The reason for this choice is, that 
the usual logarithmic broadening seems to be more frag- 
ile towards breaking of the symmetries discussed above. 
Using a Lorentzian broadening turned out to be much 
more stable. Nevertheless, NRG, being non-perturbative, 
can provide results for selfenergies and spectral functions, 
from which the ground state properties can be estimated. 
To our knowledge, the NRG has never been used for sys- 
tems with inter-site inter-orbital hopping. 

As discussed in the beginning of this section, the struc- 
tures of the selfenergy can distinguish usual topologi- 
cal insulating states from possible edge Mott insulating 
states. Therefore, we first look at the w-dependence of 
the selfenergy. In Fig. [TUl E(w,y) for U = 8 is shown as 
an example. A general feature is that the off-diagonal ele- 
ments are much smaller than the diagonal elements. Fur- 
thermore, the selfenergy for any site y can be well approx- 
imated by E(w, y) = [E(0, y) + (1 - - ib(y)oj% 

where z(y) = [1 — dT,(y)/du)] u=0 . This behavior of En ! 22 
is similar to those in correlated metals. Therefore, the 
result can be characterized as a renormalized gapped 
Fermi-liquid with renormalization factor z(y) and quasi- 
particle broadening ~ b(y)uj 2 . There are no characteris- 
tic features in E(w), such as isolated poles at |w| < Ati, 
which would drive the edge states into Mott insulating 
states. Moreover, the prefactors b(y) are significantly 
smaller than those in metals, because the density of states 
around cu ~ is suppressed by the bulk gap. This 
simple w-dependence of E(w) agrees with the previous 
DMFT+IPT calculations. DMFT+NRG confirms that 
this dependence actually holds for larger values of U than 
the band width W . The edge states always remain renor- 
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FIG. 10: (Color online) Selfenergy for U = 8 and L y = 10. 
Red, green, blue curves are En, £13,21 and E22, respectively. 



malized massive Dirac fcrmions as long as U is smaller 
than the critical value of the Mott transition U c . The 
edge Mott insulating states are not found in our cal- 
culations. Possible inter-edge Umklapp scattering does 
not seem to have significant effects on the edge states, 
although its coupling constant g^f might be large com- 
pared to the gap width A^ti when U 3> A^ti- We note 
that all correlation effects can be characterized by a single 
quantity U/W, and concerning the Mott insulator tran- 
sition, U ~ Ati does not introduce a new characteristic 
energy scale, although Ati can seemingly characterize 
the effective band width of the edge states. 

Although we cannot find the edge Mott insulating 
states in the present calculations, there is another possi- 
bility for Mott insulating states to appear. In inhomoge- 
neous strongly correlated systems, such as optical lattices 
with harmonic potentials,-- Mott insulating behaviors in 
some domains of the systems have been found at finite 
temperatures. For these states, the site dependence of 
local correlations is most important. 

In the following, we will discuss the spatial dependence 
of correlation effects by looking at the renormalization 
factors z(y). In Fig. II H we compare z(y) for several 
values of the interaction strength U. While for weak 
interactions, the y-dcpendcncc of z(y) is moderate, for 
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FIG. 11: (Color online) The renormalization factor z(y) for 
several values of U for L v — 10. 



strong interactions, e.g. U = 12 and U = 14, which 
are close to the critical value for the Mott transition, 
z(y) is largely suppressed, especially at the edge sites. 
That can be understood by recalling that the coordi- 
nation number at the edge sites is smaller than in the 
bulk, which restricts the electron motions at the edges. 
Strong correlation effects at boundaries due to the re- 
duction of the coordination number have been discussed 
in the literature for various systems, such as metallic 
surfaces , 44 ' 45 heterostructures,— ~— optical lattices^ - — 
and spin systems^ Even though the electrons are largely 
renormalizcd especially around the edges, z(y) is nonzero 
for all the sites and the bulk gap due to the SO inter- 
action remains finite. This state can be characterized 
as a renormalized TI, and therefore can be adiabatically 
connected to the non-interacting TI without closing the 
gap. In this sense, the system remains topologically un- 
changed for all interaction strengths U < U c , where U c is 
the critical interaction strength for the Mott transition. 

The orbital averaged local density of states 
— (l/27r)Im(Gn + G22) for small uj directly shows 
the renormalization effects as seen in Fig. [T2J Small dips 
inside the bulk gap correspond to the finite size gap, 
and finite values for \oj\ < Ati are due to the massive 
Dirac edge states. As has been discussed in Sec. Ill B[ 
both Ati and Af are renormalized in a similar way. 

For larger values of U, U > 14, the renormalization 
factors become zero z(y) = for all y, signaling Mott in- 
sulating states for all sites. Our results are qualitatively 
equivalent to the results for the single band Hubbard 
model where both edges and bulk are metallic without 
interactions. Although z(y) decreases especially around 
the edge sites, as has also been found in the previous 
studies 4 ^ - — , z(y) vanishes for all y at the same time at 
zero temperature, and a nearly homogeneous Mott insu- 
lating state is stabilized. This can be explained by quan- 
tum tunneling of electrons from the bulk to the edges. If 
the bulk stays metallic, then, the edges are imposed to 
be metallic even though z at the edges is strongly sup- 
pressed. 
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FIG. 12: (Color online) The local density of states at y — 1 
(red curve) and y = 5 (blue curve) for several values of U < 
U c . 



In the TIs, on the other hand, the bulk is gapful, and 
therefore, one could expect that the quantum tunneling 
of electrons is diminished, which possibly leads to inho- 
mogeneous solutions where only the edges are gapped 
with a Mott gap Ami- Such solutions could be inter- 
preted as heterostructures of Mott insulators and TIs. 
Nevertheless, the Mott insulating sites near the edges 
can reduce the effective system size of the TI, leading to 
an enhanced finite size gap. 

However, as mentioned above, inhomogeneous Mott 
insulating states are not obtained in our zero temper- 
ature calculations. The absence of such inhomogeneous 
solutions can already be understood within the linearized 
DMFT. As was shown in Figs. ElandES the local selfen- 
ergy can be approximated as £11,22(^1 y) — £11,22(0, + 
(1 — z(y)^ 1 )uj + 0(w 2 ) for small oj, and £12,21 is much 
smaller than £11,22 in magnitude and can be neglected. 
In this approximation, £(w) has the same structure as 
a simple two orbital metallic system, and the discus- 
sions in the previous studies 4 ^ -45 ' 48 ' 49 can be easily ap- 
plied to the system. Within the linearized DMFT, in- 
homogeneous solutions such as z{y = l,L y ) = and 
z(l < y < L y ) ^ are not stable provided that £12,21 
is negligible. Because the bulk gap is generated by the 
inter-site inter-orbital hopping, the gap generating mech- 
anism extends to the edge by electron motions, making 
the whole system topological insulating. In other words, 
the system is never cut by z = at specific sites when 
the interaction strength is increased. This argument is 
based upon the linearized DMFT equations and is a di- 
rect consequence of the simple structure of the local self- 
energy Z w (u>,y) ~ [£«<(0,2/) + (1 - z(y)- 1 )uj]5w which 
can be checked in the DMFT+NRG and DMFT+IPT 
calculations. We also note that, even if spatial correla- 
tions are taken into account within cluster extensions of 
the DMFT, Mott transitions for large size systems are 
expected to be homogeneous^ 

However, especially at the edge sites, the Fermi-liquid 
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FIG. 13: (Color online) The imaginary parts of the local 
Green's functions for several values of U at y = 1 (upper 
panel) and y — 5 (lower panel) when L y — 10. Red, green, 
blue curves are Gn,Gi2,2i and G22, respectively. 



behavior in our system would be fragile against tempera- 
tures, and Mott insulating domains could be seen at finite 
temperatures as it was discussed in the literature! 49 ' 50 
One can expect that there are temperature ranges where 
the peak structures around ui ~ in the local density of 
states at the edge sites disappear, while they survive in 
the bulk. In such a case, Mott insulating domains could 
grow from the edges to the bulk as U is increased. Finite 
temperature studies will be discussed in the future. 

Finally, we briefly discuss the gap structures in the 
imaginary part of the local Green's functions which cor- 
respond to the local density of states. Figure [TS] shows 
ImG for several different interaction strengths U. With 
increasing U for U < U c , the structures around the Fermi 
energy become more and more renormalized, and finally 
for U > U c a Mott gap with Ami ~ U opens. Here, we 
note that the Mott transition as a topological phase tran- 
sition is qualitatively different from transitions in which 
parameters in the non-interacting part of the Hamilto- 
nian are changed. For example, if the parameter Mq 
is changed (for U = 0), a transition from a TI to a 
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FIG. 14: (Color online) A schematic picture of the In- 
dependence of Ati and Ami- At finite temperatures or if 
spatial correlations are taken into account, Ati would not 
close and the transition would be discontinuous. 

trivial band insulator will occur. However, the gap in 
both phases is well defined in momentum-space. The 
gap width in each phase is changing continuously near 
the topological phase transition. At the Mott transition, 
however, the well-defined gap in momentum-space Ati 
for the TI phase is replaced by the Mott gap Ami which 
arises due to a gap in the atomic electron configurations 
in real space. Although the Mott transition would be 
continuous in the present calculations with the DMFT 
at zero temperature, at finite temperatures or if spatial 
correlations are taken into account, the transition can 
be expected to be discontinuous. 2 ''' 42 '' 5 The schematic 
behavior of the gap is shown in Fig. IT41 

Finally, let us make some comments on £12,21. As 
mentioned in the previous section, the local selfenergy 
satisfies £12,21 (w, y) = — £12,21 (^j y)*, and especially, 
Re£i2,2i(w = 0,y) = for all y sites. However, 
Im£i2,2i(w = 0, y) does not necessarily have to vanish 
due to symmetry arguments. Although our calculated 
£12,21 is much smaller than £11,22 and could be safely 
neglected, we have small but finite values of £12,21 = 
0,2/) near the edge sites. This might be related to the 
inter-edge Umklapp scattering. However, since £12,21 (w) 
shows complicated w-dependence within the NRG cal- 
culations, it is difficult to evaluate £12.21 = 0) with 
sufficient accuracy. 

IV. SUMMARY AND DISCUSSION 

In this article we have discussed correlation effects in 
two dimensional topological insulators. In the first part, 
we investigated finite size effects in relatively large sys- 
tems. The effective theory for those edges states is given 
by one-dimensional spin-1/2 fermions including inter- 
edge correlated forward scattering. 



The Tomonaga-Luttingcr parameters are K s > 1 for 
the spin sector and K c < 1 for the charge sector. They 
include the renormalization factor defined in the whole 
two dimensional system and might not be so much away 
from unity. The analysis of the scaling dimensions of the 
perturbations revealed that the tunneling between two 
edges is the most relevant perturbation around the non- 
interacting point when the Fermi wavenumbcr of the edge 
states is kp ~ 0. We evaluated the gap size induced by 
the tunneling for each sector which might be observed in 
experiments. 

For small systems, the edge states are massive Dirac 
fermions and there are no massless single particle excita- 
tions. We discussed the correlation effects by using the 
inhomogeneous DMFT with the IPT in the ribbon geom- 
etry. We confirmed that the naive expectation based on 
the Fermi-liquid picture of the correlated TIs, that the 
amplitude of the finite size gap is simply renormalized by 
the averaged renormalization factor, and the localization 
length of the edge states is not affected by the interaction, 
actually holds for weak interactions. 

Finally, we discussed the Mott transition using the in- 
homogeneous DMFT+NRG. The w-depcndcncc of the 
selfenergy is qualitatively similar to that in the metal- 
lic systems. We could not see characteristic features 
in the selfenergy which signal the edge Mott insulating 
states, although inter-edge Umklapp scattering is taken 
into account within DMFT. We showed that correlation 
effects are significantly stronger near the edge sites due 
to the reduction of the coordination number. Since the 
edge states are localized around the edges, the relatively 
strong correlations largely affect them especially near the 
Mott transition. As the interaction is tuned, we found 
a trivial Mott transition from the topological insulating 
state to a nearly homogeneous Mott insulating state. The 
numerical results can also be understood within the lin- 
earized DMFT as a result of the quantum tunneling of the 
electrons as has been discussed for systems without bulk 
gaps. The transition to a Mott insulator is characterized 
as a topological phase transition, and the character of the 
gap in each phase is totally different and the gap closing 
is not required at the transition in general. 

In the present DMFT calculations, we have focused 
on the paramagnetic phase and neglected instabilities to 
magnetically ordered phases^ However, in addition to 
homogeneous magnetic states, magnetic moments only 
near the edge sites might be possible. 74 ' 75 We also as- 
sumed that all parameters in the Hamiltonian are uni- 
form, although there can be changes especially near the 
edge sites in experimental situations. These issues are 
future problems. 
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Appendix A: Details of NRG 

In this Appendix, we show details of the NRG 
formalism. DMFT + NRG calculations with off- 
diagonal selfenergy have been performed in the context 
of superconductivity , 70 i 71 and our treatment is an exten- 
sion of the previous works. We also discuss details of the 
numerical calculations. 



1. Effective Impurity Anderson model 

The impurity model arising from the DMFT self- 
consistency reads 



Hiam{v) 
He. 



H c (y) + Hhyb(y) + H imp (y), 

al a (ms)gii' (T (ms)ai' <T (ms) 



ms.W <7 



H hy b = ^ \fi <T hi'(r{'ms)ap a {'i 



msAVa 



+a\ l(T (ms)h* Wa (ms)fi a , 



H», 



la 



Elcrflcrfla + ^tllll/ij, flu] 



where 



Et = (M - ft, -M - fi) 



The impurity is represented by the operators fi a , where 
I is the orbital- and a the spin index. The conduction 
electrons are given by the operators ai a (ms), where ms 
corresponds to an discretization interval (m+), (m— ) as 
defined below. Therefore, in this context gu, corresponds 
to the kinetic energy, while hw corresponds to the hy- 
bridization between impurity and conduction electrons. 
Note that all the parameters in Hiam depend on the site 
index y, which is suppressed for clarity of the notation. 
From now on, we will neglect the spin indices, because 
the lattice Hamiltonian is diagonal in spin space. The 
energy is discretized into intervals 



Im+ — (^m+1 j *£m] 



with x m = x A _Tn (m € N). 

We consider the particle-hole symmetric case and as- 
sume that gu (ms) S M. and can be written as 



g(ms) 



£(ms) C(ms) 
C(ms) — £(ms) 



(Al) 



This assumption is legitimate, because our lattice Hamil- 
tonian is defined only by real values when using open 
boundary conditions for the y-direction. The hybridiza- 
tion h is parametrized as 



h(ms) 



a(ms) b(ms) 
b' (ms) c(ms) 



where a, b, b' and c are also assumed to be real. 

The kinetic energy term g is diagonalized by the fol- 
lowing orthogonal matrices, 



U(ms) 

u 2 (ms) 

v 2 (ms) 
E(ms) 



u(ms) —v(ms) 
v(ms) u(ms) 

1 ( x + €(ms) 



1 - 



E(ms) 
£(ms) 
E(ms) 



y/P(ms)+C 2 (ms). 

E(ms) can finally be identified (for both s = ±) as 

E(ms) = \x m +a; m +i|/2, 

corresponding to the discretized intervals. The kinetic 
energy and the hybridization are transformed as 

g(ms) = g(ms)U 

E(ms) 
-E(ms] 

h(ms) = h(ms)U(ms) 

au + bv —av + bu 
b'u + cv —b'v + cu 

A B 
B' C 

Then, in this new basis a, the impurity Anderson Hamil- 
tonian is written as 



^ a) a (ms)E a (ms)a a (m,s), 



ms,a 



^2 [f}hi a {ms)a a (ms) + al(ms)hi a (ms)fi], 



ms.La 



( 



h c 

Hhyb 

a a (ms) = y^?7^(ms)a;(ms), 
l 

where E a (ms) — (E(ms), —E(ms)). 

The retarded non-interacting Green's function of the 
impurity Anderson model can be written as, 

G~ 1 (uj + iQ) =uj- E f -K(u + i0), 
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where the hybridization K is given as, and 

Kw{uj) = 2J - g{ms)]~ 1 h){ms)) w h(m+) 

ms 
m s 

\ -> \ -> hi a hl' a 

2^2^ u-E a 

ms a 

Using the identity 1 = J ded(e ~ E a (ms)) for each a, we 
obtain 

K w (w) = / de — , 

J w — e + iO 

A «'( e ) = hia{ms)h Va {ms)5{e - E a (ms)). 

We note that due to the Kramers-Kronig relation, 
A(w) = —knK(u) 

IT 

is satished. 

Next, for each interval 7™n we define 



a(m+) b(m+ 
b'(m+) c(m+ 

6(m— ) — a(m— ) 
c(m— ) — b'(m— ) 



Therefore, a possible additional condition which keeps 
the symmetry of h is, A(m—) — —B'(m+) and C(m— ) = 
+B(m+). These conditions are linear in ft.. In the limit 
C — > 0, the first equation corresponds to 6_ = — b',, 
and the second one to — b'_ = b + . Here, we have used 
shortened expressions X s = Y s > for X{ms) = Y(ms'). 
Then, the equations for h are, 





= A\ + A 2 _, 






A22+ 


= B'l + B' 2 = 


A 2 _ 


+ B' 2 , 


A12+ 


= A+B' + + A- 


B'_ 


= A_(-A + +B'_), 


An- 


= B'l+B 2 _ = 


C 2 _ 


+ Bt 


A 22 - 


= cl + c 2 _, 






Ai2_ 


= C+B+ + C. 


B_ 


= C_(C+ + B_) 



A(ms) = J du>A(u). 

Thus, for the integral region I mo + which includes only a 
single value E(ms) € {E(na)} ntT , the following equations 
hold 

A w (m+) = ^2h ll (ms)h l > 1 (ms), (A2) 

and for I mo —, 

A«/(m-) = y]hi2(ms)hi'2{ms). (A3) 



In this set of equations the upper three and lower three 
equations are decoupled. Because these equations are 
non-linear and have several possible solutions, the so- 
lution which can be smoothly connected to the limit 
C —> +0 and b, b' should be carefully chosen. The 
results are, 



S, = - 



(A 12s ) 2 



(Au a - A 22 . 
2y An s A 2 2s 



4A? 2s 



(A 



lis 



i22s J 



A 2 

^12s 



A+ 


= VA11+ + 


-5+, 


A- 




A12+ 




+ 6+ + VA22+ + $+ 


K 


= -A-, 




B'_ 


= — \/A22+ 


+ S+, 


c+ 


= -y / A 2 2- 


+ S-, 


C- 




A12- 


VA11- 


+ <5_ + V / A 2 2- + 5_ 


B+ 


= C-, 




B- 


= -VA11- 


+ 5-. 



The number of obtained equations for h given by Eqs. 
(|A"2|) and (|A"3)) is six. A u (m±), A 22 (m±), Ai 2 (to±) are 
involved and A21 gives the same equations as A12. How- 
ever, the number of independent hi a is eight and two 
extra conditions are needed for hi a to be determined. 
One possibility is that extra conditions are imposed so 
that h is connected smoothly to the trivial situation if 
we turn off the SO interaction, which is legitimate be- 
cause 1/111,22! > 1^12, 2 1 1 is expected in general physical 
situations. 

Taking the limit C — > +0 for the entries of the con- 
duction electrons of the original untransformed impurity 
model (Eq IA1[) and because £(m+) > and £(m— ) < 
can be assumed, the limiting behaviors of u, v should be 

u(m+) —> 1, v(m+) —> 0, 
u(m— ) -> 0, v(m— ) — >■ 1, of A ~ C: <5+ = <5_,j4 h 



Note that for the above solution of 8 S , 



lim lim S s = lim lim 

lis Ai2 s — >0 Ai2 s — ^0 A22s^Ai 



6 S = 



is satisfied and S s is smoothly connected to the case with- 
out the SO interaction. Because we focus on the particle- 
hole symmetric case with Au± = A 2 2zp, A 12 + = A 12 _, 
some relations can be derived from the above expressions 



-C+,B' 



-B,,A. 
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and B'_ = B_. The resulting h has the form of 
h(m+) = 



A(m+) B(m+) 
-B(m+) -A(m+) 







h(m—) 



— > 



V /An(m+) 

-v / An(m+) 

A(m-) B(m-) 
B(m—) A(m— ) 

-VAii(m-) 



VAn(m-) 







> and &, 6' 

h{ms)U^ (ms) 



where '— >' means the limit C 
0. Correspondingly, h(ms) - 
diag ( y/Auims), - ^/An(ms)) . 

We note that a condition is required for A, £>, £>' and 
C to be real quantities: 



- max{A lls ,A 22s } < 5 S < 0. (A4) 

Another condition is needed for S± to be real: 

A lls A 22s - A\ 2s > 0. (A5) 

These two conditions are essentially equivalent. The 
above derivation is based on the assumption that h 
should be connected smoothly to the one without the SO 
interaction, C, b, b' — > 0. If one assumes other smooth 
connections, such as a connection to £ — > 0, different 
expressions for A ~ C would be obtained. 



models with inter-chain hopping have never been studied 
by using NRG. The inter-chain hopping in the present 
model arises from the inter-orbital hopping, and it makes 
numerical calculations hard. 

Let us assume that this transformation is given by an 
orthogonal matrix u = (u^J, 

msa 

a a (ms) — ^ t u nmsfnvi 
nv 

Eva v a XX 
u nms u n'ms ~~ °nn'0vv', 

msa 

Eva va t: XX 

u nms u nm's' ~ °mm' S s' °aa' ■ 

ray 

Then, the hybridization term would be, 
H hy b = E [fjhi a (ms)a a (ms) + h.c] 



msla 



E fr"«( T 



f nv i + (h.c). 



We define a vector (uo)ms = Uq^ s and an inner product 
as 

"oms = hv a (ms), 



nms n ms ' 



2. Chain Hamiltonian 

In the next step we transform -ffiAM into a chain Hamil- 
tonian. We want to transform the impurity Anderson 
model to 

He ^ Envv' fnvafnv'tr 

nvv'a 

+ ^ ] tnvv' [fnvcrfn+lv'cr + f\-\-\ v i a fnvu\, 
nvv' (7 

Hhyb = E Vvv'\f-\ va fov'<j + foy'cf-lvcr], 
vv'a 

Himp ^ ] ^vf—lurjf—^^ Hjntjf— Xvai f—lvtr]- 

va 

Here, f-\ v corresponds to the former f a , 

f—lv ^ ^ figyfa • 



This chain Hamiltonian is schematically shown in Fig. 
[5] of the main text. Because the impurity part Hi mp 
is not transformed from its original expression and H c 
again includes the inter-chain hopping, we understand 
that the index v in the basis {fnva-} characterizes the 
orbitals I = 1,2. As far as we know, impurity Anderson 



and normalize u v ^ as as u£ ms = Uoms/V K«o> < )l- We > 
then, orthogonalize {uJ,Wq} keeping their symmetry to 
obtain {wq,Uq}. To fulfill symmetry requirements, the 
trigonalization must be performed always by including 
two vectors. Therefore, we assume that orthogonaliza- 
tion can be performed by 

ul new = ul + c{ulul)ul 
u 2new = ul + c(ul,ul)ul 



with one coefficient c. 



Then, 



-Anew j-.2new\ 



(ulu 2 )[l + 2c + c 2 (ulu 2 )}. If c= -0.5, (up™, aft™ ) = 
0.25 x (uq, Uq) 2 . Because and Uq have been normalized, 
|(uJ,Uq)| < 1 is satisfied. Repeating this procedure leads 
to orthogonalized unit vectors {mJ,Uq}. The hopping 

h va (ms) 



from site 



T to site is calculated from Uq° s 



as 



Vvv' ={Uq,Uq). 

Next, let us move to the kinetic term. 
H c = y^ E a (ms)a) a (ms)a a ('ms) 



EE 



E a (ms)i 



'nms n' ms 



fnvfn 
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This claims that 

^ ] ^o( ms Ki?uVms = £nuv'&nn' 
rasa 

T ^nvv' \&n' ,n+l ~(~ ^n.n' + l] 

should hold. 

II = 0, ••• ,n are obtained, we calculate e nh 
from 



0.1 



nms nms ■ 



£m>*/' = ^ E a (ms)i 
msa 

Then, we define 

^ [(£ a (ms)L^£„ w ')<„ s 



i-n—lvv 1 u n—lms ' 



and normalize each of m^'^ to obtain normalized vec- 
tors Ujj'+i- Using as initial vectors, we perform the 
orthogonalization introduced above, and obtain orthogo- 
nal vectors Usual Schmidt orthogonaliza- 
tion should also be imposed so that (Uj l+1 ,u'j ) = for 
j = 0, ■ ■ ■ , n are satisfied. Then, we evaluate t n as 



tnuv' = ^2 E a (ms)u 



VOL 

nms 



'n+lms ' 



All the parameters in the above decrease exponentially, 
and the input A is well reproduced from the chain Hamil- 
tonian derived in this way. 

We note that, although our chain impurity Anderson 
model is rather complicated, it is smoothly connected to 
the one without the inter-orbital hopping which is noth- 
ing but two decoupled single-orbital impurity Anderson 
models when there are no inter-orbital interactions at 
the impurity site. The ground state of our impurity An- 
derson model is simply the Kondo screened state or the 
decoupled state depending on the structure of the hy- 
bridization A. 



3. Discretization, Truncation and Broadening 

As already mentioned in the main text, the NRG is 
not able to gain its usual accuracy around the Fermi en- 
ergy due to the spin-orbit interaction which creates an 
hybridization between the two orbitals. For the most of 
our shown results we used a discretization A = 3.0 and 
a truncation of N K = 5000 states. Although N K = 5000 
may sound a large number of kept states, one should 
compare this number to that for a single impurity An- 
derson model (SIAM) for which similar accuracy with 
Nsiam = V5000 w 70 states would be obtained. Re- 
ducing the discretization parameter A to the often used 
value of A = 2.0 results in a breaking of the symme- 
tries discussed above, e.g. (Gi2(w) ^ G2i(oj)) which are 



kept states = 50 
kept states = 200 
A=2 kept states = 




Discretization A 

FIG. 15: (Color online) Truncation and discretization effects 
for the metal insulator transition in the ordinary one-orbital 
Hubbard model^. We show the renormalization z for U = 
1.4W w 0.94{/ c . The lines are meant as guide to the eye. 
The blue line is the reference value calculated for A = 2, 
Nk = 1000 (independent of the discretization parameter). 



supposed to be satisfied. To reduce the truncation effects 
and restore the symmetries we had to use a discretization 
parameter A = 3.0. 

Therefore we briefly want to discuss the effects of dis- 
cretization and truncations for the well studied metal in- 
sulator transition in the one-orbital Hubbard model at 
half filling^ In Fig. [15] we show the z-values z = [1 — 
dT,/du;}~^ for different truncations and discretization 
parameters A for an interaction strength U = I AW ps 
0.94U C close to the Mott transition at U c . Large dis- 
cretization parameters A lead to faster decreasing hop- 
ping parameters along the Wilson chain, reducing the 
effects of truncated high-energy states on the kept low- 
energy states, thus restoring the fundamental working 
principle of NRG. On the other hand, a large discretiza- 
tion parameter reduces the available number of states, 
and calculation results go away from those for the con- 
tinuum limit A — > 1. 

In the usual NRG iteration scheme, the renormaliza- 
tion procedure starts at high energies (large frequencies 
in the spectral functions) and iterates down to the Fermi 
energy (zero frequency in the spectral functions). Thus, 
using small discretization parameter A and small number 
of kept states, the states that would actually couple and 
change the low energy spectrum, arc truncated at high 
energies and not correctly taken into account. 

From Fig. [15j one can read off a clear tendency for the 
calculations including strong truncation effects: correla- 
tion effects are overestimated in calculations using small 
discretization parameters A, and possibily, small A even 
lead to an insulating solution. Large discretization pa- 
rameters, on the other hand, seem to underestimate the 
correlation effects. 

For checking the presented results we have performed 
some additional test calculations reducing the discrctiza- 
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FIG. 16: (Color online) Real part of the self energy for 
U — 8 between two discretization parameters in NRG and 
IPT as impurity solver. The upper panel shows y-site 1, the 
lower panel y-site 5. 



tion parameter A but increasing the number of kept 
states. Figure [T|)] shows a comparison for the real part 
of the selfenergy for (A = 3, Nk = 5000), (A = 2, 
Nk = 8000), and IPT as an impurity solver. The ten- 
dency agrees with the results for the one-orbital Hubbard 
model suffering from the large truncation effects, sec Fig. 
IT51 The correlation effects for A = 2 are much more pro- 
nounced than for A = 3. The IPT result lies between the 
two NRG results. Furthermore, increasing the interac- 
tion strength using discretization A = 2, leads to a Mott 
transition at U ps 12. However, in all the calculations us- 
ing A = 2 we find no signs of edge Mott insulating states 
and all claims made in the main text remain valid. 
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